Anderson Transition in Disordered Graphene 
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We use the regularized kernel polynomial method (RKPM) to numerically study the effect disorder on a single 
layer of graphene. This accurate numerical method enables us to study very large lattices with millions of sites, 
and hence is almost free of finite size errors. Within this approach, both weak and strong disorder regimes are 
handled on the same footing. We study the tight-binding model with on-site disorder, on the honeycomb lattice. 
We find that in the weak disorder regime, the Dirac fermions remain extended and their velocities decrease 
as the disorder strength is increased. However, if the disorder is strong enough, there will be a mobility edge 
separating localized states around the Fermi point, from the remaining extended states. This is in contrast to the 
scaling theory of localization which predicts that all states are localized in two-dimensions (2D). 

PACS numbers: 72.15.Rn, 72.20.Ee, 81.05.Uw 
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Introduction:Giaphene is the 2D allotrope of carbon in 
which carbon atoms with sp 2 hybridization are organized in a 
honeycomb lattice. Recently, the isolation of single graphene 
sheets has become possible through chemical exfoliation of 
bulk graphite til |2|, |3[], or epitaxial growth, either by chemi- 
cal vapor deposition of hydrocarbons on metal substrates J3l 
or by thermal decomposition of SiC Jsf]. The latter method 
produces graphene layers with high crystalline quality. Prop- 
erties such as, tuning the carrier types from hole to electron 
in the same sample through a gate voltage and remarkably 
high mobility of charge carriers, even at room temperature, 
due to the absence of back-scatterings, has made the graphene 
promissing from the technological point of view in building 
the carbon-based electronic devices @]. 

The electronic properties of pure graphene can be modeled 
with a simple tight-binding picture proposed by Wallace J3l 
which provides the essential band structure of graphene as a 
half-filled system with a density of states (DOS) vanishing 
linearly at the charge neutrality point. In this picture the en- 
ergy dispersion is linear in momentum near the Fermi points, 
causing the quasi-particles behave as mass-less chiral Dirac 
electrons Jit]. 

In spite of high crystalline quality, disorder can not be 
avoided in currently available samples of graphene. Various 
kinds of intrinsic disorder such as surface ripples and topolog- 
ical defects, as well as extrinsic ones like vacancies, ad-atoms, 
charge impurities affect the electronic properties of graphene. 
The presence of a minimal conductivity at zero bias yfl im- 
plies the existence of extended states at the charge neutrality 
points, in contrast to prediction of the scaling theory in 2D 

There has been various theoretical and numerical ap- 
proaches to study the effect of different types of disorder in 
graphene (as a review, see Ref. and references therein). 
Based on one-parameter scaling theory ifioll . Ostrovsky and 
coworkers found a metal-insulator transition in graphene for 
long range disorder Ml in . The hypothesis of one-parameter 
scaling was subsequently verified with numerical calculation 
of Bardarson et al. lfl2ll . They used a transfer operator method 
to investigate the localization properties of the low-energy ef- 
fective (Dirac) theory. Similar result for the beta function was 



obtained independently by Nomura et al. 11311 by evaluating 
the Kubo formula for the conductance. Although these nu- 
merical works verify the hypothesis of single-scaling theory, 
but they predict that Dirac fermions remain delocalized at ar- 
bitrary strength of disorder. The effect of roughness on the 
electronic conductivity was studied in 11411 . and it was found 
that all states remain localized in the presence of random ef- 
fective gauge fields induced by ripples. 

Lherbier et al. used a real space and order N Kubo formal- 
ism to calculate time dependent diffusion constant D(E, t) for 
the Anderson model on honeycomb lattice lfl5ll . For small val- 
ues of disorder strength W, they found D(E, t) saturates to a 
constant value in long time limit, indicating the presence of 
extended states. At larger values of W they report a decrease 
in D(E, t) pinpointing the onset of localization. 

Most of the analytical methods used in the studies of dis- 
order in graphene, are able to handle restricted regimes of 
specific types of disorder. Their predictions are valid only 
on the low-energy scales around the Dirac points, where the 
inter-valley scatterings from the impurity potential can be ig- 
nored. When disorder strength becomes comparable to the 
band-width, it is important to take into account all energy 
scales simultaneously along with possible interplay between 
different energy scales. There also remains an important ques- 
tion that, whether there is a mobility edge in graphene in the 
strong disorder regime or not? 

In this paper we use the kernel polynomial method 
(KPM) [16] based on the expansion of spectral functions in 
terms of a complete set of polynomials to accurately calcu- 
late various spectral properties, including the density of states 
(DOS). The computation time in KPM method grows linearly 
with system size. Matrix manipulations can be done on the fly, 
which reduces the memory usage enormously. Therefore one 
can study very large lattice sizes in a moderate time. Regular- 
ization of KPM method known as RKPM remedies the Gibbs 
oscillations lfl7ll . and therefore is capable to handle any type 
of disorder with arbitrary strength in essentially exact way. 
By this method, one can treat the low-energy and high-energy 
features of graphene on the same footing, and hence the in- 
terplay between the Dirac features and high energy parts of 
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the spectrum, as well as inter-valley scattering is taken into 
account. This method enables us to explore new regimes of 
disorder strength with fascinating properties. 

Model and method: From the single-particle point of view, 
disorder can affect the non-interacting electrons in graphene, 
mainly through spatial variations off on-site energy (diagonal 
disorder) or changes in the hopping integrals due to the varia- 
tions in the distances or angles of the p z orbitals (off-diagonal 
disorder). In what follows we study the graphene with on-site 
uncorrected disorder. We consider non-interacting electrons 
moving on a honeycomb lattice in the presence of local diago- 
nal disorder. The basic model to describe this kind of problem 
is the Anderson model: 
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The first term describes hopping between nearest-neighbor 
sites and e;'s in the second term are the random on-site po- 
tential uniformly distributed in the interval [—W/2,W/2], 
In Eq. (Q]i, the energy t is associated with nearest neighbor 
hopping integral, which is about 2.7 eV in graphene. This 
model has recently been studied by transfer-matrix method in 
Ref. UM < where it was found that all states are localized, in 
agreement with scaling theory of localization in 2D. However, 
our results are in contradiction with this finding. 

In our work, to investigate the localization properties, we 
employ the so called typical DOS as a quantity which deter- 
mines whether a given state with energy E is extended or lo- 
calized Hill . The return probability for an extended state at a 
given energy E is zero IU9I1 . Therefore the self-energy of an 
extended state acquires an imaginary part to account for the 
decay of return probability, while the self-energy for localized 
states remains purely real. This reflects itself in the local den- 
sity of states (LDOS): 



p s {E)=Y,\(s\E k )\ 2 8{E-E k ), 



(2) 



where \s) denotes a localized basis state at site s, while \Ek) 
is an energy eigen-vector corresponding to energy E k . Ex- 
amining LDOS at K s sites provides a measure to distinguish 
localized states from extended ones. To do so, one defines the 
typical DOS which is a geometrical average of LDOS's, 
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where K r is the number of realizations used in numerical cal- 
culations. We also need the total spectral weight which can be 
obtained from the following arithmetic averaging: 
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where D is the dimension of the Hilbert space on which the 
Hamiltonian H is acting. 



The core quantity in both Eqs. (01 and (|4]i is the LDOS. In 
order to calculate quantities of this sort, the KPM method [ 16] 
employs a complete basis set. 

The basic idea behind the KPM is to expand the spectral 
function, say, p s (E), in terms of orthogonal polynomials, 
<fi m (E), In general, all types of orthogonal polynomials can 
be used. In the case of e.g. Chebyshev polynomials one has: 
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where 



p s {E)T m (E)dE = -(s\T m (H)\s). (6) 



Here, H is obtained from H by a simple linear transformation 
to ensure that the eigenvalues of H are in [—1,1]. The same 
procedure when applied to Eq. (01 gives: 



P {E)T m {E)dE = ^Yx[T m {H)]. (7) 



Eqs. (0 and (01 can be evaluated with a recursive relation 
first discussed by Wang I20I1 . The Tr in Eq. (0 can be 
most conveniently evaluated with a simple stochastic sum- 
mation employing the recursion relation among Chebyshev 
polynomials 11611 . However, since M in Eq. (01 is finite in 
computer implementations, one faces the classic problem of 
Gibbs oscillations. There are standard attenuation factors sug- 
gested in the literature which minimize such unwanted oscil- 
lations [16, 2^1 . Due to peculiar Dirac dispersion in graphene, 
none of the so called g— factors worked. The solution around 
this problem is to use regulated (Gaussian broadened) polyno- 
mials to calculate the moments 1 1711 : 
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Our present calculation is based on using Legendre polynomi- 
als in (0 with a equal to A/M. 

Results and discussions: To investigate the Anderson tran- 
sition in graphene, we use Eq. ([1) with different values of W 
and then calculate the LDOS via the RKPM method with M 
equal to 5000 - 12000 on a lattice with 10 6 sites. The re- 
sult for different values of W is shown in Fig. As can be 
seen in Fig. Q]the average DOS for various values of W < t 
resembles that of perfect graphene. The typical DOS, pt yp , 
is non-zero every where, indicating that none of the states 
are localized in this regime. The role of disorder at weakly 
disordered regime is to slow down the Dirac quasi-particles, 
with their Dirac nature preserved. This result is in agreement 
with other works 112 ll 12211 . Fig. [2] shows the renormaliza- 
tion of the Dirac fermion's velocity. By increasing the dis- 



order strength, the slope of p decreases. Therefore, according 
to p(u) (x \u>\/vp, the vel 
when W increases II15[ Bill . 



to p(u) (x \u>\/vp, the velocity of quasi-particles decreases 
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FIG. 1 : (Color online) p (solid line) and p typ (dotted line) for dif- 
ferent values of disorder strength W in the weak disordered regime 
(W < t). In this regime the disorder only renormalizes the Fermi 
velocity of Dirac quasi-particles. Our results are stable with respect 
to change in both lattice size, as well as M. 




FIG. 2: (Color online) Average DOS as a function of energy. By 
increasing W in weakly disordered region, the velocity of Dirac 
fermions decreases. 



Increasing W to the strong disorder regime, W > t, we 
observe a mobility edge in graphene in Fig. [3] The mobility 
edge starts at the Fermi energy and keeps moving to both left 
and right side separating localized states around Fermi energy 
from those at higher energies. This is unlike the usual scenario 
of localization where states at the band edge start to localize 
first. It has to do with the geometrical nature of graphene 
lattice where states at Fermi energy correspond to momenta 
around the K points with their wave-length of the order of 
atomic separation. Upper edge of the conduction band and 
lower edge of the valance band of clean graphene arise from 
r point of the Brillouin zone which correspond to small mo- 
menta, or long wave-lengths, which naturally localize later 
than short wave-lengths modes. Such peculiar band picture 
in graphene leads to a particle-hole continuum in graphene 
which is drastically different from usual continua of Fermi 



FIG. 3: (Color online) For W > t a mobility edge starts to appear. 
Ptyp is zero around the Fermi energy. Therefore localized states fall 
into a "gap" separating extended states of the valance and conduction 
bands. 



liquids with extended Fermi surface |23']. 

The idea of band gap opening in graphene has already been 
proposed by several groups. Pereira et al 12411 by selectively 
producing vacant sites in only one sub-lattice, found the ap- 
pearance of a clean gap without any states in it. Substrate 
induced band gap, which is accompanied by breaking the 
particle-hole symmetry is another possible scenario 12511 . 

However, our result suggest the possibility of disorder in- 
duced gap in graphene. This gap is defined as the distance 
between the upper and lower mobility edges around the Fermi 
point of graphene. We find that in the intermediate disorder 
regime (W ~ t) states at the Fermi point start to localize, 
thereby opening a small gap around it (Fig. [3j. The states in- 
side the gap are localized and can not contribute in transport 
phenomena. Upon increasing disorder strength, this gap rises 
slowly up to W ~ 2t, above which its growth speeds up with 
a high slope (Fig.|4]i. This behavior can be assigned to cross- 
over from weak to strong localization regimes. Naumis also 
obtains such disorder induced gap, using a real space renor- 
malization group scheme i26ll for vacancy doped graphene. 
Nevertheless this gap continuously grows from zero as a func- 
tion of doping ratio. 

Finally with increasing W beyond W c /t = 2.5 ± 0.5, all 
states in the band become localized. The "error" ±0.5 here 
needs clarification: One might argue that the typical DOS 
is not the best quantity to distinguish localized states from 
extended ones. Our benchmark runs for the known results 
of the critical W for the 3D cubic systems gives Wj^ D jt = 
16.0 ± 0.5. Therefore the possible intrinsic error in finding 
W c in this method is less than 0.5. 

Conclusion: For the Anderson model, in weak disorder 
regime, we observe that the Dirac fermions remain delocal- 
ized up to W* w t. In this regime, the effect of disorder 
is to decrease the velocity of Dirac fermions, hence result- 
ing in a renormalized Dirac cone. However, upon increas- 
ing the disorder strength beyond W* in this model, we ob- 
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(Color online) Gap (distance between the mobility edges 
the neutrality point) versus disorder strength. 



serve a mobility edge which supports the beta function pro- 
posed in Ref. Il 111 . Our results support the idea proposed 
by Suzuura and Ando II27I1 which explains the suppression 
of weak localization in graphene. According to their argu- 
ment, the change in relative weights of the two components 
of chiral electrons wave functions induces a new Berry phase, 
when these electrons move along a closed path. In the ab- 
sence of inter-valley scatterings (weak disorder regime), this 
new phase changes the sign of the wave-function in a given 
path with respect to its time reversed counterpart, hence lead- 
ing to a destructive interference of the two paths. This is the 
reason for the existence of extended states close to the charge 
neutrality point in weak disorder regime. We think the local- 
ization of the all electronic states, obtained by transfer matrix 
method llitl . might be due to the missing of the chirality ef- 
fect. When the disorder width is comparable with the hop- 
ping energy (W ~ t), the inter-valley scatterings are present 
and the weak-localization will be eventually recovered around 
the Fermi points lEil I29I 3oL 31]. In very strong disordered 
regime (W > 2t), localization quickly spreads over all the 
energy spectrum. In this case the system is no longer homo- 
geneous in a sense that it divides into regions with different 
chemical potential and transport is described in terms of per- 
colation in the real-space B21I33I1 . 
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